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ABSTRACT 



The development of turbulent gas flows in the intra-cluster medium and in the core 
of a galaxy cluster is studied by means of adaptive mesh refinement (AMR) cosmo- 
logical simulations. A series of six runs was performed, employing identical simulation 
parameters but different criteria for triggering the mesh refinement. In particular, two 
different AMR strategies were followed, based on the regional variability of control 
variables of the flow and on the overdensity of subclumps, respectively. We show that 
both approaches, albeit with different results, are useful to get an improved resolution 
of the turbulent flow in the ICM. The vorticity is used as a diagnostic for turbu- 
lence, showing that the turbulent flow is not highly volume-filling but has a large 
area-covering factor, in agreement with previous theoretical expectations. The mea- 
sured turbulent velocity in the cluster core is larger than 200 km s~^, and the level of 
turbulent pressure contribution to the cluster hydrostatic equilibrium is increased by 
using the improved AMR criteria. 



Key virords: Hydrodynamics ~ Instabilities - Methods: numerical - 
general - Turbulence 
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1 INTRODUCTION 

I The potential importance of turbulence for the physics of 
. galaxy clusters has been widely recognised in recent years 
' (see the discussion below for details and references). The pre- 
cise nature of turbulence in the intra-cluster medium (ICM), 
■ however, remains controversial. The relevant dimensionless 
' parameter for the onset of turbulence is the Reynolds num- 
ber Re: 

LV 

i?e = ^ , (1) 

where L is the integral scale of the flow instability, V is the 
characteristic velocity at scale L, and v is the kinematic vis- 
cosity of the fluid. For a fully turbulent flow. Re > 1000. 
In the framework of galaxy cluster physics, the uncertainty 
in Re results mostly from the determination o f u. In the 
unmagnetised case, the Braginskii formulation l|Braginskiil 

[l96^ for the viscosity is generally used, lea ding to the fre- 
quen tly reported estimate of Re ~ 100 (e.g. ISunvaev et al.l 

I2OO3I ) in the ICM. 

This assumption d oes not hold in a magnetised plasma. 
As known from theory (|Spitzerlll962l ). the presence of mag- 
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netic fields suppresses the transport coefficients below the 
unmagnetised value by a factor / which depends on the 
tangled field structure. While h ighly uncertain , it is es- 
tima ted in the range 10-2-1 (|Revnolds et all l2005l : see 
also lNaravan fc Medvedevll200ll ). An effectively suppressed 
(/ ~ lO"'^) viscosity would lead to Re ^ 1 bot h in the clus - 
ter cores and in the hotter and less dense ICM (|jonesll2007l ). 
On the other hand, it has been claimed that a factor / > 0.1 
is consi stent with the observ ation of filaments in the Perseus 
cluster l|Fabian et al.|[2003lj FI. Irrespective of a more precise 
knowledge of v, it appears that the fiow in the ICM is, even 
with a conservative estimate, mildly turbulent. 

The turbulent nature of the flow in the intra-cluster 
medium (ICM) could be directly conflrmed with the help of 
high -resolution X-ray spe c troscopy of emission line broaden- 
ing l|Sunyaev et al.l l2003l : iDolag et al] |2005| ; iBriiggen etHI 



^ A weakly suppressed, or unsuppressed viscosity is also required 
for prese rving the morphologi c al stability of A GN-driven bub- 
bles (e.g. iRevnolds et al.ll2005l : ISiiacki fc Spring cl 2006) and for 
the cluster heating by viscous dissipation i FabTan e t al. 20033"^, 
though the numerical simulations of Sijacki fc Springcl (2OO3) 
question the importance of this contribution to the cluster en- 
ergy budget. 
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l2005al : iRebusco et~aLl bOOST ). The unfortunate flaw in the 
main instrument of the Suzaku satelhte postponed this test 
to the near future. Nevertheless, other observational clues 
have been interpreted as evidence for the turbulent state of 
the IC M, such as the analysis of pressure maps of the Coma 
cluster (jSchuecker et al.ll2004l ). the lack of resonant scatter- 
ing in the 6.7 keV He-lik e iron Kq line in the Perseus clus- 
ter IChurazov et al.ll2004l ). the broade ning of the iron abun - 
dance profile in the core of Perseus (iRebusco et all 120051 ) 
and other galaxy clusters (|Rebusco et all |2006|). and the 
Faraday rotation maps of t he Hydra cluster ( Vogt fc EnfilinI 
l2005l : lEnfilin fc Vogti [20061 ') . In addition to the astrophysi- 
cal problems mentioned above, an improved knowledge of 
turbulence in the ICM can shed light on the amplifica- 
tion o f magnetic fields |Dolag et al.l |2003 : iBriiggen et alj 
l2005bl : [Subramaiu an et al. 20061 ). non-thermal emission in 
clusters llBru nctti 2004 ), and acceleration of co smic rays 
l|Miniati et"al.i,2001a .bi: IBrunetti fc Lazarianll2007l ). The role 
of turbulence has also been investigated as a source of 
heating i n cooling cores, as describ ed in the analytical 
model by iDennis fc Chandrar] l|2005l ). where both the dis- 
sipation of turbulent energy and the turbulent di f fusion 
are taken into a ccount (cf. also iKim fc NaravarJ l2003l : 
iFuiita et al]|2004al lbh. In this context, the turbulence driven 
by galaxy motions in th e IC M ha s been studied e.g. by 
iBregman fc DavidI and lKimI (|2007. ). 

From a theoretical point of view, the production of tur- 
bulence in galaxy clusters has been ascribed to two main 
mechanisms. This work does not address the outfiow of ac- 
tive galactic nuclei (AGN), which infiate buoyant bubbles 
and eventually stir the ICM, but will be focused on turbu- 
lence produced during the hierarchical formation of cosmic 
structures. 

Several numerical simulations of galaxy cluster forma- 
tion show that episod es of active merging can produce tur- 
bulence in the ICM llRickeij Il998l: iNorinan fc Br-m]ll999l : 
iTakizawallioOol : iRicker fc Sarazinli200ll : iDolag et al.H2005h . 
with typical velocities of 300 - 600 km s~^ a nd injection 
scales of 300 - 500 kpc. In their simulations, iDolag et al.l 
l|2005h use a l ow-viscosity ver sion of the GADGET-2 SPH im- 
plementation llSpringel 20051) , whic h helps to better resolve 
turbulent flows. lAgertz et ahl (|2007l ) showed that grid meth- 
ods are more suitable than SPH in modelling dynamical in- 
stabilities. 

Since the properties of grid-based methods with re- 
gard to modelling turbulence are substantially better un- 
derstood than SPH, we intend to explore their capability 
in the context of cosmological simulations. These, in turn, 
rely strongly on adaptive mesh refinement (AMR) in order 
to follow the evolution of strongly clumped media. Using 
AMR for turbulent flows is a new and rapidly ev olving field 
iKritsuk et al.1 12006| . l2007l : ISchmidt etahl |2008| ) which we 
extend to cluster simulations in this work. 

The importance of a proper deflnition of the criteria 
for triggeri ng grid refinement in turbulent fiows has been 
studied by llapichino et al.l (|2008l ) (hereafter paper I) by 
testing novel AMR criteria designed for resolving turbu- 
lent fiows, and applying them in a simplified subcluster 
merger scenario. In this paper, we extend this work to cos- 
mological AMR simulations of galaxy clusters. In this frame- 
work, we investigate how the AMR criteria tested in pa- 
per I can be profitably (from a physical and computational 



point of view) used for resolving the fiow in the ICM. In 
AMR simulations of structure formation, both the refine- 
ment of turbulent fiows and the the identific ation of small 
subclumps with a relatively low overdensity (lO'Shea et al.l 
l2005bl : lHeitmann et al.ll2007^ may be equally important; the 
relative emphasis has to be determined on a case-by-case ba- 
sis. Both aspects and the related refinement strategies will 
be addressed in this work. 

We focus our analysis on one single galaxy cluster out of 
the several objects forming in the simulated computational 
volume. An extended study of turbulence should be based 
on average cluster p roperties, sampling several objects o f 
different masses fcf. lOolag et al.1 l2005l : 1vazza et al]l2006l ). 
Such a study is beyond the scope of this work and is left for 
future analysis. 

The paper is organised as follows: in Sec. [2] we present 
the setup of the cosmological simulations and the set of em- 
ployed AMR criteria. We discuss some general properties of 
the simulated cluster and present a performance comparison 
of the simulations in Sec. |3l In Sec.|4l the results are shown 
together with comparisons between the different runs. The 
results are summarised and discussed in Sec. (5] 



2 DETAILS OF THE SIMULATIONS 

Our work is based on the analysis and comparison of several 
cosmological simulations. In this section, their setup and the 
different AMR criteria are presented. The simulations were 
performed using the AMR, g rid-based hyb rid (N-Body plus 
hydrodynamical) code enzo (jO'Shea et al.ll2005a 'FI. The dif- 
ferences between the runs and the most significant modifi- 
cations to the original public source code are presented in 
Sec. [22] 



2.1 Common features 

We performed hydrodynamical simulations in a fiat ACDM 
background cosmology with ilm = 0.3, Qh = 0.04, h = 0.7, 
(78 = 0.9, and n — 1. The simulations were started with 
the same initial conditions at redshift Zin = 60, using the 
lEisenstei n fc "HuI (|i999) transfer function, and evolved to 
z — 0. Cooling physics, feedback and transport processes 
are neglected. An ideal equation of state was used for the 
gas, with 7 = 5/3. 

The simulation box had a comoving size of 
128 Mpc h-^. It was resolved with a root grid (AMR 
level Z = 0) of 128^ cells and 128^ N-Body particles. A 
static child grid {I — 1) was nested inside the root grid 
with a size of 64 Mpc 128^ cells and 128^ N-Body 

particles. The mass of each of the particles in this grid 
was 9 X 10^ Mq h~^. Inside this grid, in a volume of 
(38.4 Mpc h~'^f, grid refinement from level I = 2 to I = 7 
was enabled according to the criteria prescribed in Sec. 12.21 
The linear refinement factor A'' was set to 2, allowing an 
effective resolution of 7.8 kpc h^^ at maximum refinement 
level. 

The static and dynamically refined grids were nested 
around the place of formation of a galaxy cluster, identified 



ENZO homepage: |http://Ica.ucsd.edu/portal/software/enzo] 
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Table 1. Summary of the cosmological simulations performed for 
this work. The first column reports the name of the run, the sec- 
ond the criteria used for the grid refinement from z = 2 (discussed 
in the text) . The third column contains the number of AMR grids 
at 2 = 0. Besides the reference run A, the horizontal line divides 
the remaining simulations in two groups, as explained in the text. 



Run 


AMR criteria 




A 


OD 




2612 


B 


OD H 


h 1 


3871 


C 


OD H 


h 2 


3882 


D 


OD H 


h 1 + 2 


5358 


E 


OD, super-Lagrangian 


4100 


F 


OD, low threshold 


5340 




from 


2 = 60 





in a previous low -resolution, DM-only si mulation using the 
HOP algorithm (|Eisenstein fc Hutlll99^ . This cluster had 
a virial mass Mvir = 5.8 x 10^'' Mq and a virial radius 
Rvir = 1.35 Mpc to within 1.3% and 0.4%, respectively, 
for all simulations. At 2 = 60, the dynamically refined part 
of the computational domain contained more than 99.5% of 
the N-Body particles that were within a virial radius from 
the cluster centre at 2; = 0. 

2.2 AMR criteria and other features 

An overview of the most relevant features of the performed 
cosmological simulations is presented in Table [1] These runs 
differ with respect to the AMR criteria that were used after 
z = 2 (with the exception of run F) . Test calculations using 
the new criteria starting at 2; = 60 showed no significant 
differences. 

Until z = 2, all the simulations were run with the 
customary refinement criteria based on the overdensity of 
baryons and DM. In both criteria, a cell is refined if 

Pi > /ipo^iV' , (2) 

where po = 3H^/8itG is the critical density. In the case of 
baryons, pi — ph (baryon density) and Qi — Qh', in the DM 
case, pi = pdm (DM density) and fii = SIdm. 

The overdensity factors fi for baryons an d DM are cru- 
cial fo r the resolution of cosmic structures (jO'Shea et al.l 
l2005bl ). If they are set too high, the AMR may fail in iden- 
tifying low overdensity peaks, resulting in a deficiency of 
low-mass halos. In our simulations, unless stat ed differently, 
we set /b = /dm = 4.0. This is the same as lO'Shea et al.l 
l|2005bl ) for DM, and smaller by a factor of 2 for baryons. 
This overdensity criterion for baryons and DM is shortly 
named "OD" in Table [T] In our reference run A, only this 
criterion is used, whereas in the other simulations it is either 
modified or used in combination with other criteria. 

In Table [1] the additional simulations are subdivided 
into two groups, corresponding to two different methods for 
better refinement of the turbulent fiow. In the first group 
(runs B, C and D), we us e the AMR criteria b ased on con- 
trol variables of the flow (jSchmidt et al.|[2008l ) and already 
tested in simulations of a subcluster merger (see paper I). 
The criteria implement a regional threshold for triggering 
the refinement, based on the comparison of the cell value 



of the variable q{x,t) with the average and the standard 
deviation of q, calculated on a local grid patch: 

q{x,t)^ {q),{t)+aX^{t) , (3) 

where A; is the maximum between the average (q) and the 
standard deviation of q in the grid patch i, and a is a free 
parameter. Similar to paper I, we tested the square of the 
vorticity uj'^ {uj = S/ X v) and the rate of compression (the 
negative time derivative of the divergence d = V ■ v) as al- 
ternative or combined control variables. They are labelled as 
AMR criteria "1" and "2" in Table [TJ respectively. Prelim- 
inary tests showed that, in cosmological simulations, these 
new criteria are effective only when used together with the 
overdensity criterion. In run B, the criteria "OD" and "1", 
with threshold a = 6.5, are used. The run C has been set 
with criteria "OD" and "2", a = 6.0. In run D the re- 
finement is triggered by "OD" and both "1" and "2", with 
Qi = 7.2 and 0:2 = 6.2. 

Resolving turbulence in AMR cosmological simulations 
involves an additional requirement to tracking and refining 
turbulent fiows. Since turbulence is driven largely by clus- 
ter mergers, particular care has to be taken t o properly re- 
fine low-mass subclusters. N umerical studies (jO'Shea et al.l 
l2005bl : [Heitmann et al.ll2007l ) indicate that the smallest ha- 
los being captured by the code depend on the overdensity 
thresholds and root grid resolution, rather than on the num- 
ber of AMR levels. The second group of simulations (runs E 
and F) is devoted to explore the problem of structure res- 
olution, and to compare this approach with the new AMR 
criteria used in runs B, C and D. 

Run E has a super-Lagrangian correction to the over- 
density thresholds defined in equation ([2]), i.e.: 

Pi > /iponiiV'(^+« , (4) 

where = —0.2. The thresholds for refinement are 
low er than in run A , especially for higher AMR levels 
(cf. IWise fc Abeill2007l ) . 

The implementation of more effective AMR overdensity 
criteria from z = 2, as in run E, could produce a spurious 
suppression of small subclumps if they are not resolved from 
their formation at high redshift. In order to investigate this 
effect on the ICM turbulent motions, we performed the run 
F with the criterion "OD" and thresholds /b = /dm = 2.0, 
a factor of two smaller than run A, using these AMR criteria 
from 2 = 60. 



3 GENERAL CLUSTER PROPERTIES AND 
AMR PERFORMANCE 

The simulated cluster is widely relaxed, as can be inferred 
from its smooth spherically averaged radial velocity profile 
(Fig. [1} and from the time evolution of its mass accretion 
which indicates a steady growth due to minor mergers and 
excludes recent major merger events. In principle, this might 
be considered a poor case for studying turbulent fiows in the 
ICM because o f the lack of violent mot i ons driven by majo r 
merger shocks (|Ricker fc Sarazin|[200ll : iMathis et al.ll2005l ). 
Nevertheless, this cluster provides a useful test for the study 
of turbulent motions generated mostly by accretion of minor 
subclumps, thus isolating the role of this phase of turbulence 
production. 
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Figure 1. Spherically averaged, mass- weighted radial velocity at 
2 = in the reference run A. The cluster centre, needed for this 
analysis, was determined by the HOP tool as the location of the 
local density peak, with a maximum disagreement of 52 kpc 
(0.04 -Rvir) in the performed simulations. The analysis tool is part 
of the ENZO release. 



Figure 4. Projections similar to Fig. [3] but showing the projected 
AMR level in a zone with a side of 0.2 -Rvin centred at the cluster 
core. The colour white corresponds to the AMR level I = 7, and 
gray to / = 6. Baryon density contours are superimposed in gray. 
The six panels refer to the performed runs, as indicated by the 
letter at the upper left corners. 




Figure 2. Projection on the xy plane of a cube with a side of 
8 Mpc h~^{f^ 6 -Rvir)i centred on the simulated galaxy cluster, 
and showing the baryon density at 2 = in a logarithmic gray- 
scale. The plot refers to run A, but it is not significantly different 
for other runs. 

The morphology of the simulated cluster and its out- 
skirts is shown in Fig. (2] Several accreting subclusters and 
the related gas stripping are visible, together with some fil- 
aments (e.g. in the upper right corner). The overall struc- 
ture of the cluster does not change drastically in the per- 
formed simulations. However, interesting differences can be 
observed in the projected AMR structure (Fig. |3]|f|. The 

^ According to the AMR implementation of the ENZO code, ev- 



AMR grid follows the density (baryons and DM) distribu- 
tion in run A as well as in E and F, where the efi'ect of a 
lower AMR threshold results in a richer grid structure. As 
expected, the runs performed with the AMR criteria based 
on the regional variability of the control variables of the 
flow show that, in addition to the refinement on the density 
peaks, the grid structure is finer also in some locations which 
are correlated to the local velocity fluctuations and therefore 
to turbulent flows, rather than with the mass distribution. 
Since the analysis of the results will be largely focused on 
the cluster core (Sec. l4.3]l . Fig.|3]shows the projection of the 
AMR levels at the cluster centre at 2: = 0. This part of the 
computational domain is generally highly resolved, though 
the grid coverage changes in the different runs. 

From the point of view of grid complexity (and, conse- 
quently, of the required computational resources), our sim- 
ulations form a rather homogeneous sample. In particular, 
the AMR thresholds of the simulations B~F are tuned in a 
way to produce a number of AMR grids which do not greatly 
exceed 2 x iVgrids(^) (cf. Table [!}. 

As an indicator for the efBciency of AMR, the volume 
occupation fraction as a function of the AMR level is shown 
in Fig. [S] The covering differences between the AMR criteria 
are more significant at high refinement levels (up to one 
order of magnitude at level 7). Regarding the efficiency of 
resolutions at higher levels, the simulations can be roughly 
grouped in three sets: the very resolved runs D, E and F, 
the intermediate resolved B and C, and the least resolved 
reference run A. 



ery grid patch of the same resolution level is handled as a single 
"AMR grid" . The total number of such grids in the computational 
domain is reported in Table [Tl third column. 
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Figure 3. Projections similar to Fig. [2] but the projected AMR level is shown. The colour scale ranges from the AMR level 2 (white) 
to 7 (black like, for example, in the inner regions of the cluster). Baryon density contours are superimposed in gray. The radius of the 
white circles is r = ifvir- The six panels refer to the performed runs, as indicated by the letter at the upper left corners. 



4 PROPERTIES OF THE TURBULENT FLOW 
IN THE ICM 

4.1 Resolution of subcluster mergers 

Before focusing on the features of the turbulent flow we 
present a qualitative account of the effectiveness of the new 
AMR criteria in refining substructures in the ICM. This 
analysis is directly relevant since we are concerned with tur- 
bulence driven by inf ailing sub clumps. 

For brevity, we limit most of this comparison to the 
runs A and B. Figures |5] and [7] show some slices of the sim- 
ulated cluster and its surroundings ai z — 0. We verified 
that the results do not depend on the choice of the slice 
plane. From the velocity slice (upper right panel) one can 
easily identify an infalling subcluster with a mass of about 
2 X 10" Mq at coordinates (0.505; 0.494), located just 
outside the virial radius, and a smaller (mass of the order of 
some 10^^ Mq h-^) subclump at (0.498; 0.502) moving in 
the ICM. The projection of the motion of the first subcluster 
points roughly towards the cluster centre, while the second 
one moves to the left, parallel to the horizontal axis. 

In Fig. |6l one can observe that the refinement level in 
the region of this subclump is Z = 5. Its signature in the 
vorticity plot is not very prominent. In the divergence slice, 
only negative values, i.e. converging fiows, are shown, thus 



some of the visible structures are shocks. For example, a 
distinct feature (bow shock) is visible ahead of the larger 
subcluster, but it is less clear for the smaller one. 

The application of the new AMR criteria affects the 
smaller merger more strongly because it is not highly re- 
fined in run A. Conversely, in run B (Fig. [7]) one can see that 
the front of the subclump is refined to i = 6. The change 
is dramatic in the vorticity slice, where one can identify a 
large increase at the merger front and two parallel struc- 
tures, likely to be caused by the lateral shearing fiow. Also 
the signature in the divergence plot is clear and tracks the 
bow shock driven by the subclump. 

The larger subcluster, which is relatively well refined in 
the reference run, does not profit further from the applica- 
tion of the new AMR criteria. However, the lateral shearing 
flow is clearly visible and in run B a grid patch is added at 
(0.510, 0.496) because of the local peak in vorticity. 

It is difficult to compare the evolution of these mergers 
with simulations of moving subclusters in a simplified setup 
(paper I and references therein) . Even if the level of effective 
resolution is formally similar, the subclusters presented here 
have less details than in paper I because of their smaller size, 
and the background flow is far from homogeneous, in con- 
trast with the artificial merging setup. Nevertheless, some 
gross features of the smaller subclump resemble the simpli- 
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Figure 6. Slices of tlie xz plane for run A at z = 0, showing a detail of the simulated cluster and its surroundings. The coordinates are 
given in code units (1 = 128 Mpc h~^; -Rvir ~ 0.01). In all plots, baryon density contours are superimposed in black, (a): AMR levels. 
(b): gas velocity (corrected by the centre-of-mass velocity), in code units (200 km s~^ ~ 0.003). (c): Square of the vorticity modulus, in 
code units (cf. Sec. I4.2| |. (d): divergence of the flow, in code units. 



S 10 = 




I 



0.494 - 




Figure 8. Slice of the xz plane, at 2 = 0, for run B (same as the 
panels of Fig. (Tjl, showing the gas temperature in K, as indicated 
by the gray-scale. Baryon density contours are superimposed in 
black. 



Figure 5. Occupation fraction of the AMR levels at 2 = 0. The 
different lines refer to different simulations, according to the leg- 
end. The fraction is normalised to the whole computational do- 
main. By construction of the setup, the occupation of level 1 is 
0.125, and AMR from level 2 is allowed in a maximum volume of 
0.027. 



fied setup, like the side distribution of the vorticity. With 
the use of the new AMR criteria, it is also easier to visualise 
the mergers in the temperature slice of run B (Fig. El), where 
they show the well known cold front morphology. 

For completeness. Fig. [9] shows the slices of the AMR 
levels for the runs C-F. The position of the smaller sub- 
clump is resolved to Z = 6 only in run D, which implements 
the AMR based on the regional variability of u like B, and 
in run F, which has a lower AMR density threshold. The 
run C does not refine the subclump above I = 5, but this is 



© 2008 RAS, MNRAS 000.[HfT3l 



Simulations of galaxy clusters 7 




0.495 0.500 0.505 0.495 0.500 0.505 



Figure 7. Same as Fig. [6] but for the run B. 




0.4S5 0.500 0.505 0,495 0.500 0.505 



Figure 9. Slices of the AMR levels, on the xz plane, at 2 = 0, for the runs C-F (same as the upper left panels of Figs. |6] and O. The 
panels refer to the performed runs, as indicated by the letter at the lower left corners. 
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not necessarily an indication of an insufficient performance 
in refining turbulent fiows. More accurate diagnostics of the 
resolution on velocity fluctuations will be introduced in the 
next sections. 



4.2 Volume-filling and area-covering factors 

Using analytical models, ISubramanian et alj (|2006l ) study 
the dynamics of turbulence generation in clusters during the 
infall of minor mergers. Under reasonable estimates for the 
merger rate and assuming simple models for the geometry of 
the turbulent subcluster wakes, those authors found that the 
volume simultaneously occupied by merger wakes (expressed 
by the volume-filling factor fv as a fraction of the ICM 
volume) is relatively small, but the projections of the wakes 
on the sky plane fs cover a large fraction of the cluster area. 

Given a reasonable choice of the involved parameters, 
assuming a viscosity suppression factor / ^ 0.2 and re- 
quiring that the wake l ength is of the order of 2 7?vir, 
ISubramanian et al.l l|2006l ') report for a cluster mass of about 
10^^ Mq the indicative values 

fv ~ 0.25, fs = 0(1) 

which depend strongly on /. According to those authors, a 
combination of a high area-covering factor and a relatively 
lower volume-filling factor can reconcile the observation of 
ordered filaments in Perseus with indications for turbulence 
in the ICM of that cluster. The idea that turbulence is not 
completely volume-filling in the ICM is, moreover, a further 
motivation to use AMR in this proble m, saving computa- 
tional resources (cf. iKritsuk et aLll2006l ). 

In order to probe these features of the turbulent flow in 
the performed cosmological simulations, one needs to find a 
suitable way to track turbulence. Velocity fiuctuations are 
one of the distinctive features of turbulence, therefore quan- 
tities related to the spatial derivatives of velocity can be 
used profitably for this aim. In analogy with paper I, the 
norm of the vorticity will be used to probe the velocity 
fluctuations at the length scale allowed by the spatial resolu- 
tion, as a diagnostic for the resolved turbulent motions. The 
onset of eddy-like motions is closely related to the baroclinic 
generation of vorticity, expressed by taking the curl of both 
sides of the inviscid Euler equation: 

_.Vx(.x..)-— ^ (5) 

where p is the pressure of the gas. The second term on the 
right-hand side is nonzero if the two gradients are not par- 
allel, i.e. at curved shocks (jKang et al.l [20o3) and at the 
interface of infalling subclumps. 

This analysis was restricted to a sphere including the 
innermost 0.5 Rvii of the cluster, at z = 0. As a thresh- 
old for flagging a computational cell as "belonging to the 
turbulent flow", we chose the mass- weighted average of the 
vorticity norm in run A within the central sphere, (|aj|(y4)). 
The volume-filling factor fv is thus defined as the fraction 
of the analysis volume where \u}\ > {\ijj\{A)). If the vortic- 
ity is interpreted as the inverse of the eddy turnover time 
iKang et al.ll200'it . the chosen threshold corresponds to a 
large number (~ 50) of local eddy turnovers, assuring a con- 
servative estimate of the extent of the turbulent regions. 

The results of this analysis are summarised in Table |2l 



Table 2. Results of the vorticity analysis. The first column re- 
ports the name of the run, the second the mass-weighted aver- 
age of the vorticity norm in the sphere with radius 0.5 -Rvir- 
The maximum value of |t»;| is listed in the third column, and the 
volume-filling factor is reported in the fourth. The vorticity data 
are expressed in code units. 



Run 


(H> 


max (|a)|) 


Volume-filling 




(code units) 


(code units) 


factor 


A 


0.211 


1.42 


0.233 


B 


0.209 


4.04 


0.231 


C 


0.215 


2.05 


0.243 


D 


0.226 


2.75 


0.279 


E 


0.238 


1.59 


0.297 


F 


0.246 


1.49 


0.297 



The vorticity is reported in enzo code units. Dimension- 
ally, the vorticity is expressed as [t~^], and the time unit in 
ENZO is l/(47rGr2mPo(l + Zin)^y^^, with the meaning of the 
symbols introduced above. 

Interestingly, the volume-filling factor is not substantial 
and has a valu e somehow comparable w ith the theoretical 
expectations of ISubramanian et al. I liooe). Apart from the 
exact value of fv, which theoretically depends on many pa- 
rameters and, in our simulations, on the assumptions about 
the threshold of |u;|, the turbulent fiow does not appear very 
volume-filling {fv < 0.3 for all runs). 

When the results of the performed runs are compared, 
the features of the vorticity are quite different in the two 
groups of simulations introduced in Sec. 12.21 In particular, 
in the first group (with the exception of run D, to be dis- 
cussed) and fv are not much larger than for run A, 
but max(|ti>|) is. In the second group the opposite occurs: 
and fv are larger than run A, while max(|a;|) does 
not vary substantially. 

These results can be interpreted in terms of the dif- 
ference between runs that mostly refine on the fiow (first 
group) and those that refine on the substructures (second 
group). In the first case, especially for run B that refines 
explicitly only on vorticity, the maximum resolved value of 
|<i)| is larger, but it does not improve the overall efficiency in 
locally resolving the turbulent fiow. In the latter group, con- 
versely, the AMR does not focus on the turbulent fiow, but 
the grid captures the overdense subclumps more efficiently, 
which increase and fv by stirring the ICM. 

It is known (jPlewa fc Miilledl200lh that some spurious 
vorticity is generated at the boundary of grids of different 
refinement levels. One could therefore conjecture that the 
runs with the largest number of grids are potentially affected 
by this problem. The large fv in the run D might be partly 
due to this issue, because it has more grids than runs B and 
C. Nevertheless, we can exclude that this effect dominates 
the presented results. Evidence supporting this statement is 
provided by run E, which has A'gi-ids similar (to within 5%) 
to runs B and C, but very distinct flow properties. 

The area-covering factor resulting from these data can 
be visually inspected in Fig. 1101 for the representative cases 
of runs A and F. The turbulent flow, in agreement to 
ISubramanian et al.l l|2006l ), has fs close to unity. As ex- 
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Figure 10. Projection of the analysis sphere onto the xy plane, 
showing the area-covering factor fg. In the black zones, \u}\ > 
(|a;|(j4)) along the line of sight. The two panels refer to the runs 
A (left) and F (right), as indicated by the letter at the upper left 
corners. 



pected, the area-covering factor is larger in runs with a larger 
fv. 



Figure 11. Spherically averaged, mass-weighted radial profile of 
the rms baryon velocity at z = for the reference run A. The 
definition of the rms velocity is provided in the text. 



4.3 Turbulent features: the ICM and the cluster 
core 

The analysis of the geometry of the turbulent flow presented 
above is complementary to the knowledge of the magnitude 
of t he turbulent veloc ity. From an operational point of view 
(cf . iDolag et al.ll2005l ) , the calculation of a root mean square 
(henceforth rms) velocity implies the definition of an average 
reference velocity, in order to distinguish between the bulk 
velocity flow and the velocity fluctuations. 

In the case of spherically averaged radial profiles, the 
most natural choice is to use the average velocity {v(r)) in 
the spherical bin r ± Ar. The mass-weighted rms baryon 
velocity is then defined as 

where rrii is the mass contained in the cell i, and the sum- 
mation is performed over the cells belonging to the spherical 
shell of amplitude r ± Ar. 

The radial profile of Wrms for the reference run A is 
shown in Fig. If II The rms baryon velocity is a fraction 
of the virial velocity (Jvir = ■\^ /GMvir/2flvir ~ 960 krn . 
This is simil ar to the results of lNorman fc BrvanI (|l999l ) and 
iDolag et al.l (|2005j), though in the latter case the quantita- 
tive comparison is more difficult, because our cluster is out- 
side the mass ranges which those authors explored. 

The radial profile of Vi-ms at r > 0.1 -Rvir for the runs 
B-F does not differ significantly from Fig. [TT] This result 
could be considered a shortcoming of the adopted approach 
for the resolution of turbulent fiows in the ICM, but can be 
easily explained by the properties of turbulence and the fea- 
tures of the performed AMR simulations. In fact, the volume 
filling factor of turbulence is not very large, so any quantita- 
tive change in Vrms is likely to be washed out when averaged 
on a spherical shell. Closely related to this point, the AMR 
thresholds imposed to avoid an excessive number of grids 
particularly affect the ICM where the volume of each spher- 
ical shell is increasingly larger. 

In the small volume of the cluster core, where the level 



of refinement is relatively higher than in the ICM, the re- 
finement strategies are more effective and the use of different 
AMR criteria introduces more rel evant change s betw een the 
simulations. We notice that also in lDolag et al.l l|2005l ) (prob- 
ably for reasons related to the turbulence filling factor), the 
region where turbulent motions are most important is lo- 
calised in the inner ~ O.f i?vir. In the core region of the clus- 
ter, the analysis was not performed on radial profiles because 
this tool proved not to be robust for r < 0.07 Rvir, prob- 
ably because of sampling issues. In particular, the baryon 
rms velocity profile in the innermost part of the cluster is 
very sensitive to small displacements of the centre and to the 
calculation of (w(r)). Therefore, the mass- weighted averaged 
quantities within a sphere of radius R — 0.1 Rvir, centred at 
the cluster centre, provide a more robust diagnostic for the 
resolution of turbulent flows in the cluster core. 

The results of the core analysis are summarised in Ta- 
ble [3] For the calculation of Vrms, equation (O was applied, 
where {v{r)) is the mean velocity in the analysis sphere. The 
rms velocity is the quantity that shows the largest variation 
between the reference run A and runs B-F, with an increase 
of about 40% for run C. We notice that the highest value is 
predicted in a run belonging to the flrst group. 

Among the first group, run C seems more effective than 
run B. A possible reason is that the AMR implemented 
in run B is not designed to refine explicitly on shocks 
(cf. Paper I) , which are ubiquitous in the cluster me dium 
(jMiniati et al.ll2000l : IRvu et al.ll2003l ; iKang et al]|2007l ). De- 
spite of the larger number of grids, run D has features which 
are comparable with runs B and C, probably because of the 
larger thresholds used for its refinement criteria. 

In the second group, run F performs slightly better than 
run E. Comparing the density thresholds for refinement in 
these two runs (equations Q and ([2]), respectively, with the 
appropriate parameters) , one can see that the thresholds in 
run E are smaller for I > 5, namely in most of the cluster 
core. Therefore, the better performance of run F in the core 
is not due to a more intensive, local use of AMR, but is likely 
to be caused by the more accurate tracking of subclumps 
along the whole cluster evolution. 
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Table 3. Mass-weighted values of some hydrodynamical quantities, calculated within a sphere with R = 0.1 -Rvir centred at the cluster 
centre at ^ = 0. For the different runs (first column), the table lists the rms baryon velocity tirms (second column), the baryon density 
p (third), the temperature T (fourth), the ratio of turbulent to total pressure Pturb/^tot (fifth, defined in the text) and the entropy 
S = T/pT-\ with 7 = 5/3 (sixth). 



Run 


^rms 


P 


T 


-PturbZ-Ptot 


S (10^ X 




(km S~^) 


(10-2i5 g cm-3) 


(lO'^ K) 


(%) 


code units) 


A 


211 


1.13 


7.85 


1.37 


3.11 


B 


240 


1.05 


8.14 


1.71 


3.41 


C 


298 


1.04 


7.97 


2.69 


3.37 


D 


266 


1.07 


8.23 


2.08 


3.39 


E 


239 


1.09 


8.04 


1.72 


3.26 


F 


272 


1.05 


7.91 


2.26 


3.32 



Other variables are also listed in Table |3l For density, 
temperature and entropy the variations introduced by the 
use of different AMR criteria are smaller than in the case of 
frms (about 10%, at most, for density and entropy, and 5% 
for temperature). 

The ratio of turbulent to total pressure (fifth column in 
Table |3} is defined as: 



Ptuv 



./3 



kT/{iimp) -|-w?ms/3 



(7) 



where k is the Boltzmann constant, fj, = 0.6 is the mean 
molecular weight in a.m.u., nip is the proton mass, and T is 
the mass- weighted temperature (fourth column in Table |3]). 

The contribution of Pturb to the total pressure in the 
cluster centre is marginal, but we notice that it is increased 
(in case of run C, almost doubled) in the new runs. The in- 
crease of Pturb is consistent with the decrease in density: the 
turbulent motions introduce an additional term to the pres- 
sure equilibrium, which is balanced by a smaller baryonic 
pressure and a smaller central density. 

An increase in temperature and entropy with respect to 
run A are also observed. The mo derate increase in entropy 
is similar to what is observed in iDolag et all (|2005l ). when 
the low-viscosity version of SPH is used. In that work, the 
increase was ascribed either to a better shock modelling or 
to a more efficient gas mixing. Both features can be retrieved 
in our AMR simulations B-F, at some degree, so we agree 
with this interpretation. 

The trend in temperature is compatible with the en- 
hanced dissipation of kinetic to internal energy, although 
the increase is not closely related with the value of Vrms- We 
verified that the increase in the internal energy in the new 
runs is comparable, by order of magnitude, to the increase 
of the energy dissipation IS.E = (Aii)^r/25, where An is the 
increase of nrms, 5 is the effective spatial resolution and r is 
a time on the Gyr-scale. 

Figure [T^] displays the radial profile of the turbulent to 
total pressure outside the cluster core, for the run A. Equa- 
tion [7] was used for the computation of Pturb /Ptot(''), with 
i)rms and T averaged in the spherical shell r ± Ar. Interest- 
ingly, the profile increases with r because of the decreasing 
temperature profile and the slightly increasing velocity dis- 
persion profile. This behaviour suggests a growing impor- 
tance of the turbulent contribution to the total pressure, up 
to 10% within R = 0.5 Pvir, where the ICM can be rea- 



0.10 



0.01 L 



0.1 



0.5 



Figure 12. Spherically averaged, mass-weighted radial profile of 
the ratio of turbulent to total pressure (defined in the text), at 
2 = 0, for the reference run A. 



sonably assumed in hydrostatic equilibrium (HSE), as also 
indicated by Fig. [1] Similar to Fig. 1111 the pressure ratio 
profile is not changed by the use of new AMR criteria. Nev- 
ertheless, the magnitude of the turbulent pressure contribu- 
tion suggests that a better turbulence modelling (cf. Sec.[5| 
could have interesting outcomes in the ICM. 

Finally, we note that the DM velocity dispersion profile 
is basically unchanged along the performed simu lations, con- 
sisten t with our methods and in agreement with lDolag et al.] 
(|2005l ). 



4-3.1 Convergence test of the AMR criteria 

In order to perform a verification test on the new AMR 
criteria, it is important to check whether the increase of 
^^rms5 with respect to value computed for the reference run A 
is actually due to a better refinement of the turbulent fiow, 
or to spurious correlations with other grid-related quantities 
like A'^grids or the grid coverage of the cluster core. 

First, in the presented results we notice that there is no 
obvious correlation between the values of «rms (Table |3]) and 
the grid coverage of the cluster core (Fig. Q, meaning that 
the magnitude of the rms velocity does not depend trivially 
on the level of core resolution. Analogously, it is true that 
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Figure 13. Projections on the xy plane, showing the projected 
AMR level in a zone with a side of 0.2 -Rvin centred at the cluster 
core, at z = 0. The colour white corresponds to the AMR level 
I = 8, light gray to Z = 7 and dark gray to I = 6, respectively. 
Baryon density contours are superimposed. The two panels refer 
to the runs with maximum AMR level Zmax = 8 and the AMR 
setup indicated by the letter at the upper left corners. 



run A has less grids than the other runs, but Vrms does not 
seem directly correlated with A'grids (cf. Table [TJ . This is 
clearly shown by the comparisons of run B with C, or D 
with F, which have almost the same A^grids at z = but 
different core velocity dispersion. 

A more significant convergence test has been carried out 
by repeating the runs A and C with the same parameters 
and setup, but allowing (from z = 2) different maximum 
AMR levels. The grid coverage in the cluster core at 2 = 
is shown in Fig. [13] for the runs with maximum AMR 
level /max = 8, whereas the two runs with Zmax = 6 have a 
complete grid coverage in the selected region. 

The rms velocities within the innermost 0.1 Rvir are 
presented in Fig. 1141 For each AMR setup, one can see that 
VrmsC^max = 6) < Wrms('max — 7), because the flow structures 
are less resolved by the coarser grid structure. On the other 
hand, the runs with Zmax ~ 8 have a level of rms velocity 
which is the same as the simulations with Zmax ~ 7. This is 
also due to the AMR grid, because the volume coverage of 
Z = 8 is very limited with the employed refinement criteria 
(Fig.O. 

It appears therefore that the level of velocity fluctu- 
ations in the cluster core, when averaged on a sphere of 
R = 128 kpc , is not affected by setting the effective 
resolution below 7.8 kpc /i"^ (AMR level Z = 7). The val- 
ues of Urms show a clcar convergence for increasing Zmax. It 
is interesting to point out that iirms(^) < frms(C) even at 
/max = 6, when the grid coverage of the cluster core is iden- 
tical for the two AMR setups. From this convergence study 
we can assert that the effectiveness of the new AMR criteria 
in resolving turbulent flows does not come trivially from a 
larger grid coverage, or a larger number of grids, but instead 
from a more suitable refinement strategy. 

Though our simulations show convergence between the 
AMR levels 7 and 8, one should note that the number of 
DM particles (and thus the mass of each particle) could 
affect the results. This issue has not been fully explored in 
this work, but in a test runs with coarser particle resolution 
(64^ DM particles at Z = and 1) and AMR setup A, we 
found a velocity dispersion only 5 % smaller than in run A 
with Zmax ~ 6, which formally has the same effective spatial 
resolution. 



AMR setup A 
AMR setup C 



maximum AMR level 



Figure 14. Mass- weighted, rms baryon velocity calculated within 
a sphere with R = 0.1 -Rvir centred at the cluster centre at 2 = 0, 
for runs with different AMR maximum levels (indicated on the y 
axis) and refinement criteria (cf. Tablc[T]l. Solid line: AMR criteria 
as in run A. Dotted line: AMR criteria as in run C. 



5 SUMMARY AND CONCLUSIONS 

The study of turbulent flows in the ICM is important for 
a thorough understanding of galaxy clusters and for the 
physics of the plasma in these objects. The turbulent state 
of the ICM is, from an observational point of view, still a 
debated issue, and the theoretical determination of the kine- 
matic viscosity contains large uncertainties. Besides these 
open problems, a specific question has been addressed in 
this work: are the existing numerical techniques of grid- 
based, AMR cosmological simulations suitable for studying 
the main properties of turbulence in the ICM? In partic- 
ular, our attention was focused on the AMR criteria em- 
ployed in simulations of strongly clumped objects. A galaxy 
cluster was therefore simulated with a standard, reference 
setup, and then the simulation was repeated with five dif- 
ferent choices of refinement criteria. 

The runs were subdivided into two groups, with an un- 
derlying difference in the refinement strategy. In the first 
group, we used AMR criteria based on the regional variabil- 
ity of control variables of the fiow, designed for refining the 
grid at the locations of significant velocity fluctuations. In 
the second group, a more accurate tracking of the overdense 
subclumps was enforced. 

The distinction between the two groups of runs is im- 
portant for the interpretation of the results of Sects. f42l and 
14.31 In the runs B, C and D, the values of max(|a;|) (Ta- 
ble [2| and ?;rms (Table [S]) are larger with respect to run A, 
and to the runs of the second group. The AMR criteria used 
in these runs are therefore very suitable for refining the grid 
where the velocity fluctuations are underresolved, and result 
in an increased magnitude of the turbulence. In the second 
group of runs, the AMR does not refine explicitly on tur- 
bulence, but the lower threshold on overdense region allows 
a better tracking of subclumps. In other words, the stirring 
mechanism for turbulence generation in the ICM is better 
resolved. Therefore comparatively larger (|a;|) and a larger 
volume-filling factor are obtained. 

In both groups of simulations, the results in the inner- 
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most 0.1 -Rvir are consistent with a better modelling of the 
turbulent flow. In addition to an increased rms velocity, the 
change in Pturb /Ptot indicates that the turbulent component 
plays a non-neg l igible role in the pressure support. Recently 
IChurazov et al. (|2007l ) , based on X-ray and optical observa- 
tions of M87 and NGC 1399, put a constraint of the order of 
10 % on the non-thermal contribution to the thermal pres- 
sure, due to turbulent motions, cosmic rays and magnetic 
fields. This is not in contradiction to the findings of our 
work, because we find that the turbulent pressure is only a 
few percent of the thermal pressure at most. 

The presented comparison shows clearly that run A 
(which implements AMR criteria widely used in grid-based 
cosmological simulations) fails to reproduce both the mag- 
nitude of the rms velocity and the spatial extent of the tur- 
bulent flow. This central issue should be carefully taken into 
account in future investigations of the ICM turbulent fea- 
tures. These results suggest that computationally more de- 
manding simulations will be able to resolve more flow and 
substructures, which in turn can further contribute to the 
turbulent properties. 

Apparently the presented AMR approach to the mod- 
elling of turbulent flows in the ICM has limited utility out- 
side of the cluster core. As discussed in Sec. 14.31 this is partly 
due to the volume-filling properties of turbulence in clusters, 
and partly to the constraint on the number of produced 
AMR grids which, if too many, would make the simulation 
computationally unaffordable. 

T he presented results ar e similar to previo us simula- 
tions orman fc Brvan|[l999l : iDolag et al.l l2005l ). as far as 
the mag nitude of th e turbulent velocity is concerned. Inter- 
estingly, IkI^ (|2003) shows that the level of ICM turbulence 
driven by galaxy motions saturates around 220 km/s, which 
is not very different from our results. As stated in the Intro- 
duction, a detailed comparison of the features of turbulence 
in SPH and grid-based simulation is out of the scope of this 
paper, but we observe tha t the trends that we inferred from 
Table 13] are also present in lDolag et al.l (|2005l ). In that case, 
when the low-viscosity SPH scheme is used, the variation 
of some of the variables is more marked that in our study. 
However it should be stressed that in this work we do not 
make use of a improved numerical scheme to better resolve 
turbulence, but only of more suitable choices of AMR crite- 
ria. 

In many astrophysical problems (including galaxy clus- 
ters) the range of length scales needed to follow the tur- 
bulent cascade down to the Kolmogorov dissipation scale 
extends well beyond the grid spatial resolution limit, even 
using AMR. A strong improvement in the modelling of tur- 
bulence would be the application of Large Eddy Simulations, 
where only the largest scales are resolved and the dynamics 
at length scales smaller than t he spatial resolution is han- 
dled b y a subgrid scale model l|Sagautll200ll : ISchmidt et al.l 
l2006al lb[). Such a tool will allow to consistently evaluate the 
level of subgrid turbulent energy and will make the numer- 
ical information about the level of turbulence in the flow 
more directly available. We expect that the magnitude of 
the turbulent motions, also outside 0.1 -Rvir, could be sub- 
stantially increased. The use of this tool in the framework 
of cosmological simulations will be explored in future work. 
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